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Spin transport currents and the spin-transfer torques in vohage-biased superconducting- 
ferromagnetic nanopillars (SFNFS point contacts) are computed. We develop and implement an 
algorithm based on the Ricatti formulation of the quasiclassical theory of superconductivity to solve 
the time-dependent boundary conditions for the nonequilibrium Green's functions for spin transport 
through the ferromagnetic interfaces. A signature of the nonequilibrium torque is a component per- 
pendicular to the plane spanned by the two ferromagnetic moments. The perpendicular component is 
absent in normal-metal-ferromagnetic nanopillars (NFNFN) contacts, but is shown to have the same 
order of magnitude as the in-plane torque for non-equilibrium SFNFS contacts. The out-of-plane 
torque is due to the rotation of quasiparticle spin by the exchange fields of the ferromagnetic layers. 
In the ballistic limit the equilibrium torque is related to the spectrum of spin-polarized Andreev 
bound states, while the ac component, for small bias voltages, is determined by the nearly adiabatic 
dynamics of the Andreev bound states. The nonlinear voltage dependence of the non-equilibrium 
torque, including the subharmonic gap structure and the high-voltage asymptotics, is attributed to 
the interplay between multiple Andreev reflections, spin filtering and spin mixing. These properties 
of spin angular momentum transport may be exploited to control the state of nanomagnets. 

PACS numbers: 72.25.Mk, 74.50.-|-r 



I. INTRODUCTION 



Spin-polarized current passing through a ferroniagnet (F) 
can transfer spin angular momentum to the ferromag- 
net and exert a torque on the magnetic momenl^ii^ii^ii^. 
This provides a mechanism for manipulating the state 
of nanomagnets. For example, magnetization precession 
and reversal driven by the spin-transfer torque has been 
observed in ferromagnet-normal metal multilayers^-^. A 
promising multilayer geometry for device applications, 
such as random access memory and microwave oscilla- 
tors, is the so-called "magnetic nanopillar" consisting of 
a ferromagnet-normal-metal-ferromagnet (FNF) trilayer 
connected to normal metal (N) electrodes^i^. The typi- 
cal thickness of each layer is several nanometers, and the 
diameter of the pillar is of the order of 50-100 nm. 

Slonczevi^ski proposed a unified description of equilib- 
rium and nonequilibrium spin-transfer torque in mag- 
netic nanopillars based on the dynamics of the spin 
current^. The spin current is not conserved across the 
ferromagnetic interface, so the spin-transfer torque on 
each ferromagnetic (F) layer equals the net spin flux it 
absorbs. In thermal equilibrium a persistent spin current 
flows in the N layer between the two ferromagnets. The 
spin current is related to the formation of spin-polarized 
states due to the spin-dependent confinement potential 
of the FNF trilayer-iSiiiii^. The equilibrium current is 
responsible for the exchange torque between the two F 
layers. This torque prefers either parallel or anti-parallel 



alignment of the magnetic moments, depending on the 
thickness of the N layer^^iii. Voltage-biasing an FNF 
trilayer results in a nonequilibrium spin current in the N 
layer, as well as in the electrodes. The resulting torque 
on each F layer is proportional to the voltage. 

A scattering theory for computing spin current 
and spin-transfer torque in normal metal-ferromagnet 
hybrid structures has been developed by several 
authorsi^ii^iiiii^. Within this approach the ferromag- 
netic layers are treated as sharp interfaces between nor- 
mal conducting leads. Scattering matrices describing the 
transmission and refiection of electronic excitations, to- 
gether with their spectra and distribution functions in 
the electrodes, determine the charge and spin transport. 

Waintal and Brouwer pointed out that spin angu- 
lar momentum transfer acquires new features in mag- 
netic nanopillars connected to superconducting (S) 
electrodes^^'^°. For example, they showed that the 
nonequilibrium torque in NFNFS junctions may favor 
a perpendicular configuration of the two moments due 
to Andreev reflection at the NS interface for voltages 
below the superconducting gap, eV < Ai^. They also 
showed that the equilibrium torque in SFNFS junctions 
depends on the condensate phase difference, 0, between 
the superconducting leadsSS. These features are con- 
sequences of energy-dependent, phase-sensitive scatter- 
ing of the electronic excitations at the ferromagnetic in- 
terfaces. Fabrication and control of superconducting- 
ferromagnetic nanopillars appears to be within the ca- 
pability of current technology—. For example, the dc 
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Josephson effect in SFNFS junctions has been reported 
in Ref.^. 

In order to better understand, design and develop super- 
conducting spin-transfer devices, a quantitative theory of 
spin transport and spin angular momentum transfer un- 
der noncquilibrium conditions is needed. In this paper 
we develop the theory of noncquilibrium spin transport 
and spin-transfer torque in voltage-biased SFNFS point 
contacts, i.e. superconducting nanopillars with diameter 
smaller than the superconducting coherence length of the 
electrodes. Under voltage bias, V, the phase difference, 
(j), increases with time at a rate given by the Josephson 
frequency, toj — 2eV/h. In response, the spin current and 
spin-transfer torque oscillate with time. In general, the 
time-dependent torque excites the ferromagnet by setting 
the magnetic moments into motion, which in turn affects 
the spin transport. A full theory of the coupled dynamics 
of the ferromagnetic moments and the spin current is out- 
side the scope of this work. In this paper we assume the 
ferromagnetic moments do not change direction on the 
time scale set by the Josephson frequency, Tj = 27: /ljj. 
This assumption is justified for point contacts over a wide 
range of voltages because the characteristic frequency for 
the magnetization dynamics is set by the magnitude of 
the torque which scales linearly with the cross-sectional 
area of the contact, A. For junctions with sufficiently 
large cross section, our assumption holds only for high 
voltage bias. The quasi-static approximation eventually 
breaks down when the bias voltage is sufficiently small 
compared to the superconducting gap that the magne- 
tization dynamics is comparable to, or faster, than the 
phase dynamics. Some aspects of this limit are discussed 
in Ref.2". 

Our theoretical approach is conceptually similar to the 
scattering theory for spin-transfer torque in normal- 
metal-ferromagnet hybrids^^'^^. We develop a scat- 
tering formulation for spin transport which includes 
the effects of particle-hole coherence arising from the 
superconducting leads. Our method is based on 
the Ricatti formulation of the quasiclassical theory of 
superconductivity^i^i^i^^'^ and is explained in Sec. 
im We introduce the particle-hole coherence functions, 
the quasiparticle distribution functions and the bound- 
ary conditions obeyed by these functions at ferromag- 
netic interfaces. We also discuss the scattering matrices 
that describe the ferromagnetic interfaces; these matri- 
ces incorporate the key effects of spin-dependent scatter- 
ing: spin mixing and spin filtering. In this section we 
also summarize the method we use to solve the interface 
boundary conditions. In Sec. IIIIl we discuss the theory 
of the equilibrium torque in SFNFS contacts. We iden- 
tify the microscopic scattering processes that give rise 
to the out-of-plane torque and derive an analytical re- 
sult for the spin-transfer torque in terms of the spectrum 
of Andreev bound states at the point contact. In Sec. 
HVl we report results for the voltage, temperature, and 
misalignment angle dependences of the dc and ac com- 



ponents of the spin current and spin-transfer torque in 
voltage-biased SFNFS nanopillars. We summarize our 
results in Sec. |Vl Mathematical details of the solution 
of the time-dependent boundary conditions are presented 
in the appendix. 



II. THEORY 



A powerful formalism to study time-dependent noncqui- 
librium transport in superconducting heterostruc- 
tures is provided by the quasiclassical theory of 
superconductivity^^'^'^. This theory describes phenom- 
ena on length scales larger than the Fermi wavelength, 
A/, and on time scales long compared to the inverse 
Fermi energy, h/Ef. The central equations of the the- 
ory are a set of transport equations which govern the 
spatial and temporal variations of quasiclassical Green's 
functions. Strong localized potentials such as surfaces 
and interfaces are characterized by a scattering matrix, 
an interface "S-matrix" , which defines the probability 
amplitudes for scattering between incident and outgoing 
normal-state particle and hole quasiparticles. The inter- 
face S-matrix is the key input to the boundary conditions 
for the quasiclassical Green's function a^^'^^i^^ . 

Although these boundary conditions have been used 
to investigate the proximity effect and current-voltage 
characteristics of tunnel junctions and point contacts, 
their application has been somewhat limited. The 
boundary conditions are nonlinear, possess spurious un- 
physical solutions which have to be discarded and are 
time-consuming to implement numerically under non- 
cquilibrium conditions. An improved formulation of 
the quasiclassical boundary conditions was achieved by 
Eschrig^ for non-magnetic interfaces by reducing the 
boundary conditions for the quasiclassical Nambu- matrix 
propagators to boundary conditions for the particle-hole 
coherence functions, 7^/"* and -y^/"*, and the distribution 
functions for the particle- like, , and hole-like, , exci- 
tations (we follow the notation of Ref.'^^). These compo- 
nents are 2x2 spin matrices defined on classical trajecto- 
ries labelled by the Fermi velocity, vj; and are functions 
of spatial position, R, Fermi momentum p^, excitation 
energy e, and time, t. Note that pairs of particle-like and 
hole-like components, e.g. q and q, are related by the 
conjugation symmetry, g(p/,e) = q*{-Pf, -e)- 

The equilibrium coherence functions determine the rel- 
ative amplitude of normal-state particle and hole states 
that define the Bogoliubov exictations for the supercon- 
ducting state. The pair of coherence functions also com- 
bine to determine local quasiparticle excitation spectrum, 
i.e. A^(p/,e) = D^e [(1 + 7^7^)/(l - 7^7^)]. For non- 
equilibrium conditions, and in the clean limit, the co- 
herence functions also determine the probability ampli- 
tudes for branch conversion from electron-like to hole-like 
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quasiparticles, and vice versa. The distribution functions 
determine the occupation probabihties for the particle- 
and hole-hke excitations. We refer to both the coher- 
ence functions and distribution functions as Riccati am- 
phtudes, since all of these functions obey Riccati-type^S 
differential equations'^Q '^^i'^^i^"^ . These equations are sup- 
plemented by boundary conditions at an interface, and 
by asymptotic conditions deep in the bulk electrodes. 
For our purposes a complete set of boundary condi- 
tions for the non-equilibrium Riccati amplitudes were re- 
cently derived for spin-active interfacea^ii^^. In contrast 
to the nonlinear boundary conditions for quasiclassical 
propagators^^ the boundary conditions for Riccati am- 
plitudes are easier to solve and free of unphysical solu- 
tions. Thus, the Riccati formulation of quasiclassical the- 
ory provides an efficient method for finding the nonequi- 
librium Green's functions near spin-active interfaces^!. 

Several author s^^i^i^i^i^ have noted that caution must 
be exercised when applying the quasiclassical transport 
equations and related boundary conditions to multi- 
barrier proximity structures even for structures with di- 
mensions that are large compared to the Fermi wave- 
length. When there are multiple interfaces present the 
quasiclassical formulation can lead to errors that exceed 
the small expansion parameters that are the basis for 
the quasiclassical approximation, i.e Ij etc. 

This breakdown of quasiclassical approximation is be- 
lieved due to the formation of closed trajectories that give 
rise to constructive, quantum interference terms result- 
ing from scattering (for normal-state quasiparticles) by 
several different interfaces which invalidates the normal- 
ization condition for quasiclassical Green's functions^i^. 
However, for layered structures with two interfaces the 
quasiclassical theory yields correct results for junctions 
with many transmission channels. The quantum inter- 
ference effect, characterized by channel-dependent phase 
factors averages out when summed over many conduction 
channels. This has been verified by comparing the results 
of quasiclassical theory and microscopic Bogoliubov-de 
Gennes or Gorkov theory on the dc Josephson effect in 
double barrier SINIS junctions (I stands for an insulator 
^g^ygj.'^53j54j55_ j^ef^^ quasiclassical theory was applied 
to study nonequilibrium charge transport in diffusive 
double-barrier SINIS junctions by solving the boundary 
conditions at the two insulator interfaces. Our method 
of computing spin currents and the spin-transfer torque 
in magnetic nanopillars is analogous to that developed 
in Ref.^^. Given that magnetic nanopillars have many 
conduction channelsi^ii^, the quasiclassical approach is 
justified. 

In this article we apply the Riccati method to study time- 
dependent transport in superconducting-ferromagnetic 
nanopillars (SFNFS junctions). The junction geometry 
and coordinate system are shown in Fig. [1] Two single- 
domain, nano-scale ferromagnets, F° and F**, are sepa- 
rated by a normal-metal spacer, N. The FNF trilayer is 
connected in a circuit with electrodes (S) made of con- 




FIG. 1: Schematic of superconducting magnetic nanopillar. 
The left electrode, middle spacer, and right electrode are also 
called region 1, 2, and 3, respectively. The two ferromagnetic 
layers are labeled by indices a and b. Their magnetization 
directions, fia and fib, are in general mis-aligned; fia is along 
the z axis, while fit is at polar angle i/' in the yz plane. The x 
axis is along fia x fit, and the diameter of the pillar is d ^ ^o. 



ventional s-wave superconductors. The ferromagnetic in- 
terfaces, F'' and F^, are characterized by their scattering 
matrices, S° and S^, respectively, that couple normal- 
state excitations in adjacent electrodes. 

For simplicity we are assume the ferromagnets to be 
identical except that their magnetic moments are mis- 
aligned by an angle ip. The thickness of each layer, as 
well as the diameter of each nano-pillar, is assumed to be 
much smaller than the superconducting coherence length, 
^0 — fivf /2TTkBTc, and the inelastic mean free path, 
These criteria are easily met, for example, with Al elec- 
trodes. The N layer thickness L is assumed to be larger 
than the Fermi wavelength but still thin compared to the 
superconducting coherence length. A/ <C i ^ ^o- We 
also assume the leads and the trilayer are in the clean 
limit and the scattering at F" and F'' is specular. We 
choose the magnetization direction, fia, as the quantiza- 
tion axis {z axis) for spin. The magnetic moments, fca 
and /ifc, span the yz plane, while the x axis is perpendic- 
ular to both fLa and fib- As discussed above, we assume 
fia and fLb are static on the Josephson time scale, Tj. 

Figured] shows the Riccati amplitudes defined on a set of 
scattering trajectories near the ferromagnetic interfaces. 
The Riccati amplitudes in the superconducting electrodes 
are classified into "incoming" Riccati amplitudes denoted 
by lower case symbols, {'^f''^, if'^i •'^f '^f}' ^^'^ "outgo- 
ing" Riccati amplitudes denoted by upper case symbols, 
{rf ^,ff ^,Xf ,lf }, j = 1,324-34. The amplitudes in 
the N region are also labeled in this way with "incom- 
ing" and "outgoing" groups defined with respect to the 
interface with F°. 

Proximity effects, i.e. the appearance of superconducting 
correlations in the N layer, as well as the suppression of 
superconductivity in the electrodes by the FNF trilayer, 
are encoded in the local particle-hole coherence functions. 
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FIG. 2; Coupled scattering trajectories (directed lines) and 
the corresponding Riccati amplitudes for an SFNFS structure 
in the clean limit. The arrows indicate the direction of the 
Fermi momentum for quasiparticles defined on each trajec- 
tory. The right (left) going trajectories in the N region are 
labeled by 2> (2<). 



Charge and spin currents in each region, i = 1,2,3, 
are determined by the local distribution functions to- 
gether with the coherence functions. For SFNFS junc- 
tions under a constant voltage bias the phase difference 
across the junction evolves according to the Josephson- 
Anderson relation, 4>{t) = ojj = 2eV/h. Thus, all lo- 
cal Riccati amplitudes defined above are time-dependent. 
However, we can neglect spatial variations of the Ric- 
cati amplitudes within the N region. For example, the 
coherence function, r|^, in the N layer varies in space 
as r^(z') ~ r^(0)e''=^^', where kr = 2e/vf < kj is 
the Tomasch wave vector and z' is the distance from F° 
along trajectory 2>. Since z' is typically of order the N 
layer thickness (several nanometers) the phase accumula- 
tion krz' is negligible, except for grazing trajectories that 
are irrelevant for most transport phenomena across the 
trilayer. Similarly the distribution function, , varies 
as X^{z') - (O)e'"^'/''^ where ui is of the order of 
the Josephson frequency toj <^ Ef, and barely changes 
except for grazing trajectories. Thus, for thin N layers 
determining the local Riccati amplitudes is reduced to a 
simultaneous solution to the boundary conditions at F° 
and F''. 



Boundary conditions 

Boundary conditions at F" and F'' connect the "outgo- 
ing" Ricatti amplitudes to the "incoming" Ricatti ampli- 
tudes via the interface scattering matrices^. For voltage- 
biased junctions these boundary conditions describe the 
spin-dependent inelastic reflection and transmission of 
quasiparticles, loosely referred to as "multiple Andreev 
reflection" (MAR). The summation of the multiple scat- 
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TABLE I; Amplitudes for reflection and transmission of elec- 
trons and holes for the processes contributing to Eq. [T] h (e) 
stands for hole (electron), and > (<) represents right (left) 
moving. For example, \h >—> le < means a right moving 
hole in region 1 converts into a left moving electron in region 
1. 

tering processes are contained in the boundary condi- 
tions, which for the retarded Ricatti amplitudes are 

Ff =r?o7f(^?,)t+t^o7«(^?2)t, (1) 
rf =r2-o7«(5^2)t+tgo7f(^^Jt, (2) 

72^' ^r'2or^{Sl,y+4o^^iSl,y, (3) 

=4oj,"isl,V + 4or§{sl,)^. (4) 

Matrix multiplication in spin space is implied, and the o 
operator denotes the folding product, which in the time 
domain is defined as, 

/ + 00 
dt A{ti,t)B{t,t2)- (5) 
-oo 

These boundary conditions can be expressed in terms of a 
set of rules for computing the branch-conversion scatter- 
ing amplitudes near the interface. Table |T] lists the scat- 
tering amplitudes that determine the outgoing Ricatti 
amplitude, Ff , for electrode 1. This is the amplitude for 
retro-reflection by the interface with branch conversion 
from a hole-like to a particle-like excitation. The dia- 
grammatic representation of the lowest order elementary 
scattering processes that contribute to Ff are shown in 
Fig. El 

The transmission (t) and reflection (r) amplitudes in 
Eqs. [TJll are calculated from the normal-state, branch- 
conserving scattering amplitudes {S^j^.S^-^) and ele- 
mentary branch-conversion amplitudes (7^, 7^^) by sum- 
ming the amplitudes for all multiple scattering processes. 
In particular, for transmission and reflection from the 
central electrode (2) through F"/'', 

t2 - [{821)^ ~ Piio{Pl,)-\S^^^)^Y' , (6) 
= -t^2 0Piio{Pi2r\ (7) 

4 - [is'2,)^ ~ Pi, o iPt2r\s',,y] , (8) 

r'2 = -4°/333°(/332)"\ (9) 
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one lead depends on elementary scattering processes that 
couple excitations between leads 1 and 3. For example, 
for lead 1 involves the summation of multiple scatter- 
ing processes from both interfaces. Low-order processes 
are shown in Fig. |4] which include transmission and re- 
flection from both interfaces. Note that the processes 
that couple lead 1 to 3 are contained in Eqs. [T][51 

Boundary conditions for the distribution functions are 
also expressed in terms of the dressed reflection and 
transmission amplitudes, 





— 

— '1 


xf- orl - 


f- xf 


-ra 

ti 


-al 


oif oa-,(12) 




— „a 

— '2 


X2 °r2- 


F xf 


-ra 

ot.2 




oif oag, (13) 




- r'' 

— '2 






t2 


-a\ 


o^f oa^,(14) 




- 


x^ ° ^3 ^^ 


- 1\ 




-4 


olf oa^,(15) 



FIG. 3: Diagrams for interface scattering with branch con- 
version from hole to particle. The full scattering amplitude 
includes Andreev reflection and transmission with branch 
conversion to all orders in the normal state transmission 
(<S'i2,S'2i) and reflection amplitudes (S'11,5'22). 



where the coefficients /3f ■ are defined by 



Qa 



„R 



1,2. 



(10) 
(11) 



The coefficients for interface F'', /?! 



2,3, are defined 



similarly. Notice that r|*' and f |^ are the "incoming" 
amplitudes with respect to F**, so 'y^ and 7^ in Eq. [TUl 
[TT] should be replaced by F^ and F^ in the definition of 
. Expressions for t J and r j are obtained from Eq. [5][7] 
by index exchange 1 <-!■ 2, while i§ and rg are obtained 
from Eq. [8]l9] by index exchange 2 3. The boundary 
conditions for the remaining coherence functions, such as 
rf and 'J2', are straightforward to write down by evoking 
particle-hole symmetry and are not listed here^^. 



However, the distribution functions also depend on the 
Andreev transmission amplitudes, i.e. transmission with 
branch conversion, 



{r^m2-Shi^)oip-2)-\ 
(rf^2i-^fi7f)°(/3fi)-\ 



— (r3'.^32 ^32^2 ) ° {P22) ^- 



(16) 
(17) 
(18) 
(19) 



They are related to t and r by a set algebraic identities^i. 
For example, the amplitude a° describing the event 
2h <^ le < is related to by = S^j^ - aj7^5'fi. 
Note that both "retarded" amplitudes, r, t, and the corre- 
sponding "advanced" counterpart, denoted by r,t, enter 
the boundary conditions for the distribution functions. 
These advanced functions are related to their retarded 
counterparts by symmetry, i.e. q{e, e') = q{e', e)^. 



Solution to the boundary conditions 




pa pi ' p:< pi 



FIG. 4: Elementary scattering processes contributing to Fi 
originating from transmission and reflection from both inter- 
faces of the FNF trilayer. 

The left and right side conducting leads are coupled 
through the FNF trilayer. Thus, Ricatti amplitude for 



The SFNFS nanopillars are small compared to the super- 
conducting coherence length. This means that the sup- 
pression of the superconducting order parameter by the 
exchange fields of the ferromagnets, as well as the cur- 
rent flow, scales as \CA/(,a and can be neglected. Also the 
voltage drop occurs essentially at the contact because of 
the large Sharvin resistance. Thus, the values of the "in- 
coming" Riccati amplitudes in the exterior superconduct- 
ing leads, {7^, 7^, xj^, ij^; j = 1,3}, are determined by 
the spectrum and distribution functions for bulk quasi- 
particle excitations in the left (1) and right (3) leads. 
Our task is then to solve the boundary conditions, Eqs. 
[TlfT9l for the Riccati amplitudes in the N layer and the 
"outgoing" Riccati amplitudes in the exterior leads. 

We start with the coupled Eqs. and solve them by an 
iterative method. Typically we are interested in the spin 
current in junctions with moderate to high transmission 



6 



probabilities, since otherwise the spin-transfer torque ef- 
fect would be too small to be of interest. Therefore, a 
reasonable initial guess for and 7^ is obtained by ne- 
glecting normal reflections at F" and F''. Thus, to zeroth 
order in 5'"/^ Eqs. [13] yield 



(20) 
(21) 



An improved approximation for F^ and is obtained by 
substituting the zeroth-order solution into the right side 
of Eqs. [mSl After the coherence functions in the N layer 
are calculated, the outgoing coherence functions in the 
leads, Ff and F|^, are calculated from Eqs. [T]and[4l The 
transmission and reflection amplitudes (r, t, a) can then 
be evaluated from Eqs. [51fTT]and [TBlfTOl Each iteration of 
this procedure incorporates higher-order scattering pro- 
cesses. Operationally the iteration is stopped when the 
results for F|^ and 7^ converge within a specified preci- 
sion. 

The linear equations for the distribution functions, Xi^ 
and X2 , Eqs. [T^lfHl are easily decoupled. For example, 
obeys the equation 



with c — r2 o 



X^ ~ CO X^ oc - 
2 and 



d, 



o X 



3 ° ^2 - a§ o 2:3 



(22) 



(23) 



-(2^*)C' 



x^ 



~X^ 



where C-^'^ = diag[(l-7«'^7«'^)-i, -(l-7«>^7«'^)-i]. 
In particular, the spin current in the N layer is given by 



2 desin2e [ de 



Tr\^cr[g^2 



K 

> 



5f<] 



(26) 

where the trajectories are indicated in Fig. [1] and Q = 
arccos(v/ • m) is the angle between the Fermi velocity 
and the interface normal m at F", a = {(Tx, cry, az) are 
Pauli spin matrices, and Tr{. . .} denotes a trace over spin 
states. 



Interface S matrix 



The key inputs to our theory are the scattering matri- 
ces of the ferromagnetic layers, S° and S^. In principle, 
these matrices can be calculated from ab initio micro- 
scopic theory^. However, the more efhcient procedure 
is to fit a few key transport measurements to theoreti- 
cal predictions based on models for §"/''. Thus, we pro- 
ceed by carrying out calculations with basic models for 
the S''/^ matrices which capture the key features of spin- 
active scattering at ferromagnetic interfaces, i.e. spin 
filtering and spin mixing. More detailed discussions on 
the interface S matrix can be found in Ref.^''. 



Eq. [22] is solved iteratively to yield 



= 1 , 



For junctions of moderate to high transparency c is small, 
and the series summation converges rapidly. The other 
distribution functions in the N layer, X^ , , and X2 , 
can be evaluated in a similar way. The distribution func- 
tions in the leads, Xf^ and X^ , are then calculated from 
Eq. [12] and [H 

The main technical difficulty in implementing the itera- 
tion algorithm described above is the evaluation of the 
folding products. Following Ref4^ we work in the en- 
ergy domain. The folding products are then represented 
by products of infinite dimensional matrices for multi- 
ple Andreev reflection (MAR). We truncate these MAR 
matrices to evaluate the folding products. Details of our 
analysis are provided in the appendix. 

The observable spin current is calculated from the qua- 
siclassical Keldysh Green's functions. These functions 
can be constructed from the solutions for the Riccati 
amplitudes'^'*. 
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(25) 



We parameterize in terms of spin-dependent trans- 
parencies, Z?| and -D|, and a spin mixing angle, iJ^. In 



(24) particular, for interface i^". 



Qa Qa 

Oil '^12 

qa qa 

'-'21 '^22 



e-'"t'Sr iSt 
iSt e"f'Sr 



(27) 



The factors, e with (f> ~ kfL/{pf ■ m), describe the 
phase of normal-state quasiparticles propagating through 
the N layer"^^, while the normal-state spin-dependent re- 
flection and transmission matrices Sr^t are given by 











VD^e'-^/^ 

\/:D|e~'^/2 



(28) 
(29) 



The reflection and transmission probabilities for the 
two spin cigenstates are related by, Ra = 1 — Da for 
a —^,1-, defined with respect to the polarization axis 
fia of interface F". The nonequilibrium spin-transfer 
torque in normal-state magnetic nanopillars is predomi- 
nantly determined by the spin filtering effect, i.e. ^ 
we show below, spin mixing becomes 
important for the nonequilibrium torque in supercon- 
ducting magnetic nanopillars. The spin mixing angle. 
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measures the relative phase shift between spin up and 
spin down electrons upon transmission or reflection; is 
easily shown to be the angle of rotation of the electron 
spin polarization around the magnetization direction af- 
ter the electron is transmitted or reflected. This process 
is analogous to Faraday rotation in optics, and in the 
case of superconducting leads contributes to the torque 
on the ferromagnet. Note that for simplicity we assume 
the spin mixing angle for reflection and transmission are 
the same. This is a feature of models of the F layer which 
possess inversion and time reversal symmetry (including 
the inversion of the magnetic moment), see Refi^. 

Since we assumed that F^ and F° are equivalent except 
for the orientation of their moments we can obtain the 
S-matrix for contact F**, from simply by rotating the 
ferromagnetic moment by angle tp with respect to the x 
axis (see Fig. [T|), 



of, 



where U{4>) 



^'^[/(V)^.^/(^)^ si 

iU{tP)StU{iP)^ 



-2i4> < 



.6 

'22 J 



^32 



(30) 
(31) 



-i<T^if>/2 



is the appropriate spin rotation 



operator. We also note that the scattering matrices for 
holes are also related to the scattering matrices for elec- 



trons by symmetry relations^l, 5y^^''(p||) = [S'"/''(-p|| )]"■", 
where the superscript, T, denotes the transpose of a spin 
matrix and P|| = p/ — (p/ • m)m. 

In general the values of {D-^, -Dj, i?} depend on the quasi- 
particle trajectory (p/), the specific ferromagnetic ma- 
terial and the thickness of the F layers. For example, 
estimates for the transparencies of the Cobalt layer in 
Co-Cu-Co nanopillars are ~ 0.7 and Di - 0.2^^. We 
are most interested in spin-current transfer in relatively 
high transmission junctions, and thus relatively weak fer- 
romagnetic materials such as permalloy^"^^. We consider 
two representative sets of spin-dependent transparencies: 



I. Z?! = 0.81, Di = 0.64, 
II. Z?! = 0.95, Di = 0.60. 



(32) 
(33) 



We allow the spin mixing angle to take values rang- 
ing between and 2tt, with ■& — n corresponding to 
the strongest spin mixing. This is intended to model 
F layers of various thickness. To put these numbers 
into perspective, consider the delta function scattering 
model as an example of weak spin mixing. In this model, 
each F layer is approximated by spin-dependent delta 
function potentials of strength Va for spin a. In terms 
of 9a = arctan[ft(v/ • m)/VQ], the transparencies are 
Da = sin^ 9a , and the spin- mixing angle is given by 
-d = 9i — 9^ — arcsin a^/D| — arcsin y^D^. Thus, for 
transparency model I, the delta function barriers yield 
= O.OeiTT (11.0°); while for model II, d = O.U6tt 
(26.3°). 



In the S matrices defined Eqs. [27ll3T| spin flip scatter- 
ing is absent since the ferromagnetic layers are treated as 



static macrospins (single domain ferromagnets with ho- 
mogeneous magnetization) which give rise only to elas- 
tic scattering of normal-state electrons and holes. Our 
model is in line with that of Refsji^ and^S. However, our 
formalism can be generalized to extend beyond the static 
macrospin approximation and to take into account in- 
elastic scattering and spin flip processes. However, these 
considerations are outside the scope of this report. 



Transparent spin-mixing 



In the special case where normal reflection and spin filter- 
ing at F" and F** are negligible, the boundary conditions 
are greatly simplified and can be solved analytically. The 
nearly transparent spin- active point contact, introduced 
in Ref.— , highlights the effect of spin mixing on spin cur- 
rent transport in magnetic nanopillars. As discussed in 
Ref.^^, this is a reasonable approximation if the F layers 
are made of metallic ferromagnets with exchange energy 
h <§; Ef. With this approximation, the normal state 
electron scattering matrix of F° is described by only one 
parameter, the spin-mixing angle. 





gilT,1?/2 







(34) 



and the scattering matrix for F^ is given by S*" — 
Uitpp^UHi;). 

For transparent interfaces, the boundary conditions for 
Riccati amplitudes simplify, and we easily obtain the co- 
herence functions in the N layer. 



rf = g*'T.'3/2^Kg-i<Tzi?/2^ pfl ^ g-j<T,^/2~Kgi<T,^/2 
R _ icr^i>/2 R -icr'i)/2 ~R _ -i<y'-»/2~R ia^i)/2 



72 



73 e 



72 



73 e 



,(35) 



where cr^ = Uz cos-0 — ay sin-0 and cr^ = U-^. Similarly 
we obtain the distribution functions in the N layer. 



■iCT^i5/2 



'z-S/2 



= pi<yi,i)/2^K^-ia^-d/2 ~K ^ ^-ia'^f) /2 ^ia'^f) /2 



J -^2 



(36) 



III. EQUILIBRIUM TORQUE 



The equilibrium torque in superconducting nanopillars 
with chaotic N spacers has been studied in detail in 
Ref.— . Here we make explicit connections between the 
equilibrium torque and the Andreev bound states in bal- 
listic junctions. The objective is to reveal the microscopic 
scattering processes responsible for the equilibrium spin 
current, and to facilitate our subsequent discussion on 
nonequilibrium torque. 
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Analytical expressions for the equilibrium spin current 
for the N layer of the nearly transparent SFNFS nanopil- 
lar are obtained from the Riccati amplitudes. In particu- 
lar, the equilibrium Keldysh Green's functions are simply 
related to the thermal distribution function and the re- 
tarded Green's functions for the trajectories 2< and 2>, 
and g2<,^ so the spin current 



92> 



I - NfVfA- 



— tanh — Im/Tr[CT(g^ 



R 



92<)]), (37) 



where Nf is the density of states at the Fermi energy, w/ is 

the Fermi velocity and (...) = /J^^^ sin 6 cos 8(. . .)dQ/2 
is the angular average over the Fermi surface, including 
the projection of the quasiparticle spin current along the 
interface normal, i.e. m ■ Vf = cosQ. The retarded 
Green's functions are easily constructed from the Riccati 
amplitudes given in Eqs. [35l By inspecting the poles of 
g2<^ , we observe that a pair of Andreev bound states form 
at sub-gap energies 



sgn 



A. 



(38) 



with X = arccos [l — 2 sin^ cos^(?/'/2)] . For trajectory 
2>, the bound states are at e^{(f>) — e^(— </>). Thus, 
for a set of coupled scattering trajectories there are in 
general four bound states, {ef , j = 1,2,3,4}. These 
Andreev bound states are spin-polarized (c.f. Ref4S for 
a discussion of spin-polarized Andreev bound states in 
SFIFS junctions). For example, the bound states with 
energies correspond to spin eigenstates polarized along 
the axis n = (rix ^ny^riz), 

, 1p "0 4' \ / \ 

n — (sin — sin v, — sin — cos cos — cos v) . (39) 

The Andreev bound states can be visualized as polar- 
ized electrons (holes), confined within the FNF trilayer 
by the superconducting pair potential of the two elec- 
trodes. Quasiclassically, these bound states correspond 
to excitations traversing closed trajectories between the 
superconducting leads. Consider an electron on trajec- 
tory 2< for example. In order to form a bound state, 
its energy has to be such that the total phase accumu- 
lation during a round trip, including the phase acquired 
while going through the F layers as well as that acquired 
during the Andreev reflection, is a multiple of 2t: (the 
Bohr-Sommerfield quantization condition); additionally 
its spin has to be polarized along n so that in one round 
trip its spin polarization returns to its original direction. 
The evolution of the spin polarization during one round 
trip is illustrated in Fig. [5] It follows from the geometry 
that the mutual torque between F" and F** is along the 
X direction. 

In equilibrium the spin current only flows between the 
two F layers; it vanishes in the superconducting elec- 
trodes because the total spin of a Cooper pair is zero 




FIG. 5: The round trip of the spin polarization of a left- 
moving electron in the N layer with energy e!j: and initially 
polarized along n. Upon traversing F° to the left supercon- 
ductor, s is rotated from n by — 19 around fia- After the An- 
dreev reflection at the left superconductor, the retro-reflected 
hole passes through F° again so that s is rotated by another 
—I? around fia and arrives at n'. As a result, spin angular 
momentum fa = (n — n')h/2 is transferred to F", which as 
shown is along the x direction. The right-moving hole passes 
through and is retro-reflected into an electron by the right 
superconductor. After the electron passes through F*" to get 
back in N, s is rotated by — 2i? around fib and returns to n. 
The second half of the trip transfers spin angular momentum 

Tt = -Ta to F^ 



(here we only consider conventional s-wave, spin-singlet 
superconductors). The net contribution from contimmm 
states to the spin current in Eq. [37]is zero, so the equilib- 
rium spin current is carried entirely by the spin-polarized 
Andreev bound states. We carry out the integration over 
the spectrum using the residue theorem and find that the 
y and z components of the net spin current vanish, and 
that the x component is directly related to the spectrum 
of bound states. 



(40) 



where /(e) = l/(e'^/'^'^'^ -f- 1) is the Fermi function. Writ- 
ten explicitly, the spin current in the N layer is given 
by 



I{(j),ip,T) =i (NfVfATrhA) sinip x 
sin^ ■& 



-— > sm — 

2smx„f^, 2 



tanh — ■ cos 

2T 2 



q;0 



(41) 



Its magnitude is set by the scale NfVfAAh, where A is 
the energy gap at temperature T. The dependence of 
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the spin current on tjj = arccos(/ia • fib) is generahy non- 
sinusoidal since x is also a function of ip. The mutual 
torque between the F layers, 7% = — Tq = / = tx, is along 
the X axis, and if the moments were free to rotate, these 
torques would drive a mutual precession of fca and jli,. 

A similar analysis for the equilibrium Josephson (charge) 
current, 1^, can be carried out with the result that le 
may be expressed in terms of the phase dispersion of the 
Andreev bound state energies. 



It follows from Eqs. SO] and SH that 
dip d(f) ' 



(42) 



(43) 



This relation between the equilibrium spin-transfer 
torque and the Josephson charge current was previously 
written down in Reffi^ using free energy arguments. 
While this relation provides a nice check on our equi- 
librium results, it does not generalize to non-equilibrium 
charge and spin transport. 



Equations [38ll4T] are the central results for the equi- 
librium spin current. We have explicitly shown that 
besides carrying Josephson supercurrent the Andreev 
bound states mediate a mutual torque in the x direc- 
tion between the F layers. As shown in Fig. [5l finite spin 
mixing at the two F layers, i.e. rotation of the spin polar- 
ization by the exchange field, is crucial for producing the 
equilibrium spin-transfer torque. The analysis outlined 
above can be carried out for arbitrary SFNFS junctions 
with finite spin filtering. The details of the bound state 
spectrum will be different, but Eqs. I38II41I remain valid. 

The equilibrium spin-transfer torque in SFNFS junctions 
and NFNFN junctions share a similar origin. In both 
cases, the persistent spin current flow can be traced to 
the formation of spin-polarized bound states and under- 
stood as a quantum interference phenomena. However, 
there is a major difference. In NFNFN structures, the 
spin current is rather fragile to mesoscopic fluctuations^^, 
and strongly suppressed by random impurity scattering. 
In reality, the exchange coupling between the F layers in 
a normal junction is effectively short-ranged, for exam- 
ple, it vanishes for L > 20 nm^. By contrast, the spin 
current in high transmission superconducting junctions is 
robust against impurity scattering^S, and the proximity- 
induced exchange coupling is long-ranged. The long 
range nature of the induced exchange field in supercon- 
ductors by a proximity contact with a ferromagnet was 
studied in Ref,44. As demonstrated in normal metal- 
superconductor hybrid structures'*^''^^, the phase coher- 
ence between incoming electrons and Andreev reflected 
holes persists throughout the diffusive normal metal for 
energies less than the Thouless energy. These low energy 



excitations in the N layer can mediate exchange coupling 
between the F layers even when L is comparable to the 
superconducting coherence length. 



IV. NONEQUILIBRIUM TORQUE 



An applied voltage across an SFNFS junction leads to 
time-dependent non-equilibrium spin currents in the su- 
perconducting leads, as well as the interior normal-metal 
layer, i.e. currents in all regions of the trilayer, /'^''(t), 
i G {1,2,3}. Since the Josephson phase is varying as 
ip{t) = 2eVt/h = LOjt we can expand the currents in 
each layer in Fourier series defined by the Josephson fre- 
quency. 



(44) 



^ (f^}{V) cosikujjt) + f^}{V) sinikoj jt) 



fc=i 



It follows that the spin-transfer torque on F*", T^{t) = 
/(2)(t) - /(3)(i), has the form 

OC 

^{t) - ^ + E (^c coaikoj jt) + r*, ainikoj jt)) . (45) 

k=l 

For the coordinate system shown in Fig. [1] the dc torque 
^0 — '''ox^ + ■^Oyy + ■^Oz^ can be decomposed into the in- 
plane (II) and out of plane (_L) components. 



J). 

Tn 



royY - 



(46) 
(47) 



Note that both f^'" and f^'"'" are perpendicular to jib- In 
what follows we discuss how the Fourier components of 
Tb depend on the voltage bias the spin mixing angle 
and the misalignment angle ip. 



D.C. spin current 



Fig. [S] shows the dc spin current for the N layer in 
the zero-temperature limit for the delta function scat- 
tering model with = 0.81, Di = 0.64, i9 = O.Oeivr, 
and Ip = tt/2. Note that the spin current and spin- 
transfer torque are plotted in units of NfVfAA{T)h/2 
in all figures. We note two major differences between 
the nonequilibrium dc spin current in superconducting 
versus normal-state magnetic nanopillars. First of all, 
in normal-state nanopillars the spin current is linear in 
the bias voltage^- (within the quasiclassical approxima- 
tion), whereas for superconducting nanopillars the spin- 
current- voltage characteristics are highly nonlinear, and 
possess subharmonic gap structure (SGS). Secondly, the 
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FIG. 6: The x, y, and z component of dc spin current in the 
N layer for perpendicular configuration ij) = n/2. D| — 0.81, 
Di = 0.64, i9 = O.OeiTT, T = 0. The spin mixing is rather 
weak, and the subharmonic gap structures are close to voltage 



Vn 


= 2A/ne. The unit for spin current is NfVfAAh/2. 
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FIG. 7: The x component of dc spin current in the N layer 
for ip = 7r/2 and different spin mixing angles. D| = 0.81, 
Di = 0.64, r = 0. 

spin current in normal-state nanopillars is predominately 
due to spin filtering and lies within the plane spanned by 
/ta and /If, (see Ref^iS). By contrast, in superconduct- 
ing nanopillars, spin mixing is important; it controls the 
current-voltage characteristics of the in-plane spin cur- 
rent, and leads to a finite out-of-plane spin current along 
X = /ia X /if,. As shown in Fig. [6l the out-of-plane spin 
current component, Iqx, is finite at V = 0, and decays to 
zero at high voltages, eV ^ A. On the other hand, the 
in-plane spin current components, I^y and Iqz, are zero 



at y = and grow linearly with V at high voltages. 

The SGS of the spin current develops somewhat analo- 
gously to the SGS of dc charge current in voltage-biased 
Josephson ]n-sictions^^'^L^>^^. The SGS structure 
reflects the opening of additional channels ( "onsets" ) as 
well as resonances associated with multiple- Andreev re- 
flection (MAR) processes. The n-th order MAR channel 
has an onset voltage Vn = 2A{T)/ne. This transport 
mechanism is more efficient for charge transport since 
the Cooper pairs carry charge 2e, but zero spin. A sin- 
gle MAR process may transport multiple units of charge, 
but only one quantum of electron spin. Nevertheless, the 
sub-gap spin current develops due to MAR, and can be 
resonantly enhanced when a MAR trajectory in energy 
space intersects one of the FS surface bound states, the 
energy of which depends on the spin mixing angle^^^. 
This generates resonant transmission of the spin for volt- 
ages tuned to the bound state energy. For weak spin mix- 
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FIG. 8; Schematic diagram illustrating the trajectory of a 
sub-gap excitation undergoing multiple Andreev reflection 
(MAR). Charge transport occurs at each Andreev reflection, 
while spin is transported by the excitation with energy above 
the gap that escapes into the superconducting lead. 



ing, shown in Fig. [6l the SGS exhibits singularities in the 
differential conductance (d/o /dV) located approximately 
at Vn- Increasing the degree of spin mixing, increases the 
magnitude of the current, eventually leads to reversal in 
the current direction and generally smooths out the sin- 
gularities. The spin current vanishes as the spin- mixing 
angle approaches tt. Fig. [7] shows the out-of-plane dc 
spin current in the N layer as a function of voltage for a 
range of spin mixing angles. For such high transmission 
junctions, the main contribution to Iqx comes from con- 
secutive spin rotations when electrons (holes) undergo 
multiple Andreev reflection. Roughly speaking, the total 
number of subgap MARs is inversely proportional to V] 
as a result I^x decays rapidly at high voltages. The mag- 
nitude of lax also becomes small for d close to or tt, 
where the quasiparticle spin direction is hardly changed 
in each Andreev reflection (see Fig. [5]). 
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FIG. 9: The dc torque on F*' in x (upper panel) and z (lower 
panel) direction for t/; = 7r/2 and different spin mixing angles. 
= 0.95, Di = 0.6, T = O.STc 



D.C. spin-transfer torque 

Once the spin current is known in each region, the torque 
on F*" is obtained by f* = /(^J - Fig. [5] shows the 
dc torque on F** as a function of the bias voltage for 
= 0.95, Dx = 0.6, V = 7^/2, T = O.ST^, and a va- 
riety of spin mixing angles. Note that the y component 
of f* vanishes for ip — tt/2. We observe that the volt- 
age dependence of the in-plane torque (tq^) and the out- 
of-plane torque {tq^) are sensitive to spin mixing, and 
more importantly, exhibit very different voltage charac- 
teristics. The in-plane torque vanishes for V —^ 0, while 
the out-of-plane torque is finite at zero voltage. The out- 
of-plane torque varies more dramatically with the bias 
voltage for V < 2A/e, changing sign around 2 A in the 
case of weak spin mixing (i.e. i? < 7r/3). This marks 
the crossover from the low-voltage regime, where only 
high-order MAR processes with finite number of sub- 
gap Andreev reflections contribute to the spin momen- 
tum transfer, to the high-voltage regime, where direct 
transmission without Andreev reflection dominates the 
spin transport. Indeed, for V ^ 2A, both components 
become linear in V, however, note that the magnitude 
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FIG. 10: The dc torque on F** as a function of misalignment 
angle tp for strong spin mixing, i9 = 27r/3 and tt. Z)| = 
0.81, Di = 0.64, T = 0.5Tc, V = 1.5A/e. The x, y, and 
z component of the dc torque are plotted separately. Also 
shown is the magnitude of the in-plane torque (labeled by 
"in"), the vector sum of the y and z components. 



of Tg^ (for V ^ 2A) is suppressed from its normal-state 
value because of the superconducting gap in the exci- 
tation spectrum. Finally, we note that that if we scale 
the torque in units of NfVfAA{T)h/2 and the voltage in 
units of A(T)/e then both in- and out-of-plane compo- 
nents of the torque are nearly temperature independent. 

The spin-transfer toque vanishes in the parallel and an- 
tiparallel configurations for the two F layers, i.e. V = 
and TT. The functional forms for t^' "(■(/') and f^'''~(V') are 
generally complicated. Two calculations are shown in 
Fig. [TO] for d = 27r/3 and = tt at voltage V = 1.5A/e, 
with D| 0.81, Di = 0.64, and T = O.ST^. In the case 
of d — 27r/3, the magnitude of out-of-plane torque is 
considerably larger than the in-plane torque, and it pos- 
sesses a pronounced maximum around -0 ~ 7r/5. In the 
case of = TT, however, the out-of-plane torque is of the 
same order as the in-plane torque, and it varies rather 
slowly with ?/'• At higher voltage, eV > 2A, and weak 
spin mixing both the in- and out-of-plane torque vary ap- 
proximately as sinip^ as shown in Fig. [11] for = 0.06l7r 
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FIG. 11: The dc torque on F*" as a function of misalignment 
angle ij^ for weak spin mixing, i9 = O.OGItt, at high voltages, 
V = 2.5A/e. D| = 0.81, Dy = 0.64, T = O.ST^. 



and V = 2.5A/e. The high-voltage asymptotic behavior 
of 7^'" (V') is similar to that of the normal-state magnetic 
nanopillars discussed in RefiiS,. 



A.C. spin-transfer torque 



Figure [12] shows the first Fourier components of the 
ac torque on F'', fic{V) and fis{V), for = 0.95, 
Di = 0.6, § = O.UGtt, V = 7r/2, and T = O.ST^. The ac 
spin current comes from the interference between MAR 
processes of different order. As a result, the ac torque 
rapidly decays to zero at high voltages. To get a bet- 
ter understanding of the ac torque, in what follows we 
derive its analytical expression for transparent SFNFS 
point contacts. 

The nonequilibrium spin current in transparent point 
contacts can be computed directly from Eqs. [35ll36l The 
analysis is carried out in the appendix; the resulting 
expression for is given by Eq. IA17I We observe 
that the dc spin current vanishes in all three regions, 
« = 1, 2, 3. This is a feature of the nearly transparent 
interface model, and it shows that finite spin filtering is 
required for a non-zero dc torque. For the ac current the 
dominant contribution comes from the first Fourier com- 
ponent, {V) cos{u;jt) + (V) sm{u}jt). The cosine 
term is negligible compared to the sine term when the 
inelastic rate is small; we have assumed Fin = 10~'*A. 
The sine component of the current vanishes in the su- 
perconducting leads, however I^^ is finite and polarized 



along the x. Thus, the time-dependent torque is 



(48) 




eV/A '° 



FIG. 12: The first Fourier component of the ac torque on F*. 
Shown are Tie and ris in x and z direction (the y components 
are zero for Tp = tv/2). D| = 0.95, D| = 0.6, i9 = 0.1467r 
(weak spin mixing), -0 — it/2, and T — 0.5Tc. 




From Eq. IA17[ we find the torque is proportional to 



FIG. 13: Dimensionless factor, F{V, T), for the torque in Eq. 
l49lfor r = and T = O.STc F{V, T) is inversely proportional 
to V for eV > A. 



sin d and suiip, 

Tis{y,T)^NfVfAA{T)^F{V,T)l^sin^'t)) sinV'. (49) 

The voltage dependence of t* is given by the dimension- 
less function F{V,T), plotted in Fig. [fslfor T = 0.5Tc 
and T — Q. Since F{V, T) depends weakly on T, the tem- 
perature dependence of t\, is mainly determined by the 
gap A(r). 

Note that the ac torque decays monotonically with volt- 
age in transparent contacts. This is in sharp contrast 
with the ac torque, tI^{V), shown in Fig. [T^lfor < 1- 
In this case resonant MAR leads to sharp structure in 
rfgX for V < 2A. At high voltage, ey>A, the domi- 



13 




0.4 , 0.6 
tO^ t /27I 



FIG. 14: The time evolution of the spin-transfer torque on F 
in the adiabatic limit. ip = it/2, T = 0. Although tis sin{ujjt) 
is the dominant term in the Fourier expansion of T^t), the 
deviation from the sm{ujjt) dependence is obvious. 



Eq. [35] indicates that the occupation of bound states at 
energy e^{(j){t)) is given by /(A) while the occupation of 
bound states at energy €^{(l){t)) is given by /(—A). This 
procedure yields 

^W-^^/«MA^|tanhA 



sin ■& 
sinx 



|sin(| +cjt)| 



^t)\ 



(52) 



The time evolution of the torque, T^(f), at zero temper- 
ature is shown in Fig. [14] for ip = 7r/2, i? = 7r/2 and 
7r/4. If we expand Tb(t) in Fourier series, we find the sine 
part of the first Fourier component is Tis = ^NfVfAAh 
for Ip — -d — tt/2, which agrees with our numerical re- 
sult in Fig. [13] It is clear that the higher order Fourier 
components are nonzero, but their magnitudes are much 
smaller than tis. 



nant contribution to spin current comes from the lowest 
order MAR processes, i.e. the terms m = 0, 1 in Eq. 
IA17I Keeping only the leading order terms we obtain 
the asymptotic behavior for F{V, T) for Da = 1 given by 



F{V, T) ~ 2Re j 



^°°|(^f(ei)73^(^-i)7f(^) 

-xf(e)7f(e2)7i'(ei)). (50) 



At zero temperature integration of Eq. [SD] yields 

Ti^(T = 0) ~ NfVfAl\h^^(^sm^d^ sin?/', (51) 

and thus, the ac spin-transfer torque is inversely propor- 
tional to the voltage bias in the limit ey^A. 

In the opposite limit, ey<cA, the time evolution of spin- 
transfer torque is governed by the nearly adiabatic dy- 
namics of the Andreev bound states. The spin-transfer 
torque, however, does not assume its instantaneous equi- 
librium value, i.e. f{t) ^ Tcq(0(O)- The reason is that 
the dynamics of the bound state spectrum, i.e. e^{(j){t)), 
leads to out of equilibrium population of the Andreev 
bound states. The occupation of the Andreev bound 
states remains constant as a function of time until the 
bound state energy evolves to the gap edge at which point 
these states equilibrate with the continuum quasiparticles 
in the leads within a time scale much shorter than the 
period of a Josephson oscillation. This argument was 
made by Averin et al. to describe the adiabatic time 
evolution of the charge transport in high transmission 
Josephson point contacts'^^i^. To obtain the adiabatic 
limit for the dynamics of the spin-transfer torque, we re- 
place the Fermi function /(cj) in Eq. [JD] with /(±A). 
For eV > the Josephson phase ip increases with time. 
The phase dispersion of Andreev bound states given by 



SUMMARY 



In this paper we studied spin momentum transfer in su- 
perconducting nanopillars in the clean limit using mod- 
els for the interface scattering matrices for the magnetic 
contacts. Our method presented in section [II] however, 
is quite general. Further exploration of the parameter 
space {Da,'d,ip} as well as investigations on the effects 
of impurity scattering in the N spacer, spin-flip scatter- 
ing at the ferromagnetic interfaces, etc., can be carried 
out within the current theoretical framework. The avail- 
able information on the scattering parameters for SFNFS 
structures, e.g. the directional dependence of spin mixing 
angles, is a potential limitation. However, spin transport 
is dominated by trajectories close to normal incidence. 
For trajectories substantially away from the normal di- 
rection, e.g. the grazing trajectories Q ~ 7r/2, the effec- 
tive thickness of the F layer is large and the transmission 
probability is very small. Their contribution to the total 
transport current is negligible even though the spin mix- 
ing might be strong. Therefore we believe that the basic 
features we obtain from our S-matrix models for the F 
interfaces, such as the out-of-plane torque, its nonlinear 
voltage dependence, etc. will survive the average over 
scattering trajectories for more detailed models for the 
interface S matrices. Indeed we expect our model to pro- 
vide a good approximation for realistic device behavior. 

In summary, we investigated the nonequilibrium spin- 
transfer torque in voltage-biased superconducting mag- 
netic nanopillars. Our work extends earlier research on 
equilibrium phase-sensitive spin-transfer torque by Wain- 
tal and Brouwer— to nonequilibrium junctions. Our re- 
sults, and the general theoretical framework, provide a 
computational formalism for understanding the opera- 
tion of non-equilibrium superconducting spin-transfer de- 
vices. We have shown that a superconducting nanopillar 
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is an interesting system with new physics, particularly 
a rich dynamics resulting from the interplay of multiple 
Andreev reflection, spin mixing, spin filtering, spectral 
dynamics of the interface states and the Josephson phase 
dynamics. The nonequilibrium torque at finite voltages is 
in principle observable by monitoring the magnetization 
dynamics, as described in Ref.^5,. However, the details 
of the coupled magnetization dynamics and the time- 
dependent spin-transfer torque is a topic for a future 
study. 



Acknowledgments 



We thank Dr. Tomas Lofwander for stimulating discus- 
sions on this work, and for the diagram notation of the 
interface boundary conditions used in Sec. |TT1 



APPENDIX A: THE FOLDING PRODUCT 



Our evaluation of folding products follows closely the pi- 
oneering work of Arnold"*^. Double time correlation func- 
tions can be Fourier transformed between the time do- 
main and the energy domain, 



°° de' 



de" 



J-oo J-oo 

(Al) 

The folding product is most conveniently evaluated in 
the energy domain. 



[AoB]{e\e") 



(k 
27r 



A{e',e)Bie,e"). (A2) 



We can fix the condensate phase and electrical potential 
of the left electrode to be zero. It is easy to show that 
the matrix elements of "incoming" Riccati amplitudes in 
the electrodes are given by 

7f(e',e) - 2^5(e' - 6)7f (6) , (A3) 

xf(e',e) - 2n5{e'-e)x^{e), (A4) 

73^(e',e) = 2^<5(e'-6 + 2c.)7,f(e-c.), (A5) 

a;f(e',e) = 27r(5(e' - e)a;f (e + tj) . (A6) 

where uj = ujj/2, and 'j^{e) etc. take their bulk equi- 
librium value^i. The remaining Riccati amplitudes are 
obtained from those above by the symmetry relation, 
q{pf,e\e) = q* {—pf,—e' ,—e). Note that 7;f and are 
ladder operators in energy space, while the amplitudes 
Ii^tIi 'j and are diagonal in energy space. 

Inspection of the boundary conditions, Eqs. [T1ITO1 shows 
that all the coherence functions in the N layer, as well as 
all the effective scattering amplitudes (r, t, a) defined at 



F" and F**, have the form of 

+00 

A(e',e)= 27r(5(e' - e^)A„o(e) , (A7) 

m— — 00 

where we introduced the short-hand notation, = 
e -|- miu. This expansion refiects the fact that quasiparti- 
cles of energy e are only scattered to side-band energies, 
{em}, during multiple Andreev reflection processes. Fur- 
thermore, the matrix elements of the distribution func- 
tions in the N layer are nonzero only between states cor- 
responding to the side band energies, 



/(i€. 
— 2TrS{€' -€2m,)2n6{e" -e2k)X2m,2k 

(A8) 

The integration variable e provides a natural reference 
energy to introduce a discrete basis, {I'm) = \em)}, in 



energy spacer**. 



1, A„ 



[m\A\n] 



(A9) 



In the Arnold basis, each spectral operator is represented 
by an infinite dimensional matrix, and the folding prod- 
uct reduces to matrix product. 



[AoB]mn{<^) = ^A„,fe(e)Bfe„(e). 



(AlO) 



For general interface scattering matrices, the iteration al- 
gorithm described in Sec. [Til to solve the boundary con- 
ditions must be carried out numerically. To this end, we 
truncate the infinite dimensional matrices in energy space 
and work in a space of dimension (2M-I-1) x (2M-I-1). This 
corresponds to neglecting contributions from MAR pro- 
cesses with more than M reflections^^. Usually M>4 
yields sufficiently accurate results except for small volt- 
ages, e.g. eV < O.IA. 

Once the Riccati amplitudes are calculated using the it- 
eration algorithm, we can assemble the Keldysh Green's 
fmictions to obtain the spin current, e.g. in the N layer. 



M \k+m\<M 



A2mut 



I^^)it,V)^NfVfA7:-mc J2 E 

k=-M m 

X I ^(Tr[aK(fc.fc+,„)(6)] ),(A11) 



where the kernel is given 

K(m,„) (e) EE ^ ( 2m [5^ - | 2n ) 



(A12) 



The spin currents in the left and right electrodes, I*^^^ {t) 
and I^^'^{t), are computed in similar fashions. 

For transparent SFNFS point contacts, the spin cur- 
rent can be obtained analytically. In constructing the 
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Keldysh Green's functions, the matrix elements of (1 — 
A)~^ are computed by expanding it in geometric series, 
(1 — A)~^ = J^'kLo ■ ^'-"^ example, the probability am- 
plitude of bulk quasiparticles being scattered from energy 
e to e + 2fcw via MAR is 

M,(e) = (2fc|(l-rf 72«ri|O)= (A13) 

We adopt the convention nj=i(---) = 1 11^=1 = 
AkAk-i---A2Ai. Besides M^, the following transition 
amplitudes are required to construct and 

Nfc(e) = Mfe(e)e^-^''/27f(e)e— (A14) 

k 

i=i 

X7f(e-2,+2)e^^^''/^ (A15) 
Pfc(e) = Ofc(e)e^^^''/273^(ei)e-*'^-*''/2. (A16) 



The magnitudes of these amplitudes decay rapidly with 
increasing MAR order, i.e. k, as well as the bias voltage. 
The spin current in the N layer is found to be 

^ +00 oo 

m= — 00 k=0 

xy'^(Sp[?K(fc,™)(e)]), (A17) 

K(fe,™) = N,(e)i3^(6_i)N^^„(e)-Mfe(e)a:f(e)M^^„(6) 
+ Ofe(6)x3^(6i)OL„(6) - P.(6)if (e2)PL,Je) . 

/^^^ and l'^^^ have the same form as Eq. IA171 but with 
slightly different definitions of the transition amplitudes 
{Mfe,Nfc,Ofc,Pfc}. 
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